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We apply a new technique for ab initio phase determination [Acta Crystallogr. Sect. A, A55, 
48 (1999)] to solve for the average structure of the icosahedral (i) phase of AlPdMn. After an 
introduction to the crystallographic phase problem and a description of the method, we present a 
brief report of our findings for the structure of i(AlPdMn). Despite the use of data from extremely 
high quality samples, we find strong evidence of disorder in the structure, lending support to the 
random tiling model of quasicrystal stabilization. 



o 
o 
o 

(N 

c 

In 

On 



> 
On 

O 

o 
o 

c3 



o 
o 



61.44.Br, 61.10-i, 61.43.Bn 

While X-ray diffraction is one of the most useful tools 
we have for solid-state structure determination for every- 
thing from simple metals to protein crystals, in a very real 
sense diffraction experiments only provide half of the in- 
formation we need. X-ray diffractometers measure only 
the amplitudes of the various Fourier components of the 
electron density, not their relative phases. The recovery 
of the phases, given only the amplitudes, is termed the 
"phase problem", or "phase retrieval", in a more gen- 
eral context. Efforts to solve the phase problem can be 
broadly grouped into those that use additional data and 
those that do not. In the first category are three-beam 
interference measurements jj] and techniques such as iso- 
morphous substitution and anomalous scattering at mul- 
tiple wavelengths Methods of the second type are 
termed "direct" methods and were pioneered by Haupt- 
man and Karle fl. Algorithms which implement direct 
methods are highly evolved and very successful in solving 
structures of up to a few hundred atoms |5j. Recently, 
a technique called "Shake 'n Bake" (SnB) [|J|] was de- 
veloped expressly for application to more complex struc- 
tures such as proteins. Common to all direct method 
algorithms currently in use is an iteration scheme which 
modifies phases to reproduce, cither the statistical dis- 
tribution of phase invariants (triplets, etc.) expected of 
an ideal gas of equal point scatterers, or known general 
features of the electron density (e.g. "solvent flatten- 
ing"). Due to their iterative nature, the success rate of 
these algorithms, starting from random phases, is very 
unpredictable and in the case of complex structures fre- 
quently only a small fraction of a percent Q. One pos- 
sible reason for this behavior is that like all iterative op- 
timization schemes || these methods can become stuck 
in fixed points that are locally but not globally optimal. 
This phenomenon is known as "stagnation" in the image 
processing literature [ |To[ . 

Unfortunately, current direct methods are completely 
unsuitable for solving i(AlPdMn) and all other quasicrys- 
talline structures for two reasons. On the one hand, 
the properties of the quasicrystalline state require that 



their Bragg peaks be indexed by more than three in- 
dices; icosahedral phases, in particular, have an average 
structure that requires six dimensions (6D) to represent 
periodically . Furthermore, an atomistic three dimen- 
sional electron density is obtained by cutting a density 
in the higher dimensional space. The latter need not be 
concentrated at points; rather, one expects the analogues 
of atoms — "atomic surfaces" — to be rather extended 
in the additional dimensions. The inability of established 
direct methods to deal with these features of quasicrys- 
talline charge distributions recently led one of us jl2| to 
introduce a new method for phase retrieval, called the 
principle of minimum charge (minQ). This new method 
differs from previous direct methods in three important 
respects: it is not iterative or stochastic, it does not rely 
on atomicity of the charge density, and it is applicable 
in arbitrary dimension. A direct method that works in 
arbitrary dimension is useful not only for quasicrystals 
but also in image reconstruction. 

While an extended discussion of the foundations of 
minQ is available elsewhere , the central idea can be 
stated very simply: the correct phases are such that the 
total charge Q is a minimum. This claim has its roots in 
the following observation. The minimum average charge 
density (equivalently, Q) required to restore positivity 
to a tentative density reconstruction is determined by 
the deepest minimum in the electron density; the large 
number of free parameters (phases) implies that the mini- 
mum Q is achieved when a correspondingly large number 
of minima in the electron density are exactly degenerate 
and represents the relatively large fraction of the unit cell 
volume — in a correct reconstruction — having negligible 
charge density. 

Past methods to solve the phase problem in quasicrys- 
tals have been indirect and have relied on observations 
about the diffracted intensities [[l3| or relationships of the 
actual structure to closely related approximant structures 
There has also been some use of the maximum- 
entropy method Jl5[ in quasicrystal structure determi- 
nation IM, but this work has focused on structure re- 
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finement and not ab initio phase recovery. Our selection 
of i(AlPdMn) as a proving ground for the minQ method 
had a threefold motivation: as previously stated, existing 
direct methods are not applicable, important structural 
details were heretofore unknown, and a wealth of diffrac- 
tion data for extremely high quality samples is available 

0- 

The Fourier series for the electron density, for crystals 
as well as quasicrystals, may be written as 

p(r) = F + 2 cos(q • r - <£ q ) (1) 

qGA+ 

where A + is a set of nonzero reciprocal lattice vectors 
corresponding to the x-ray data (each element represent- 
ing a {q, — q} pair). F is the average charge density 
whose minimization is the goal of any algorithm based 
on minQ. We assume the size of A + is large enough so 
as to make truncation errors negligible. If we special- 
ize to the case of structures with centrosymmetry, as in 
i(AlPdMn), p(r) = p(— r) and <fiq can be either or ir. 
We may therefore replace the phases in the argument of 
the cosine with overall multiplicative signs s q = ±1. Our 
separation of Fo from the rest of the diffraction ampli- 
tudes in Eqn. (|l]) emphasizes its practical difference from 
other amplitudes in that F is not measured. Given a 
maximum spatial frequency in A + , we can evaluate the 
density on a finite grid Q. Consequently, we only deal 
with a matrix of numbers F with entries 

F rq = 2F q cos(q ■ r) reG q S A+. (2) 

Imposing positivity at our set of grid points, the min- 
imum charge principle immediately leads to a linear op- 
timization problem: 

minimize : Fo (3) 
subject to : F + F • s > (4) 
Sq e{-l,l} (5) 

This type of problem is referred to in the optimization 
literature as a mixed-integer program (MIP) |lq , ^9| : the 
objective function and constraints are linear in the un- 
knowns, while some of the unknowns are required to take 
discrete (integer) values. This general type of problem 
has a long history in operations research, with applica- 
tions to scheduling, job assignment, and set partitioning 
to name a few examples [ pof . The traveling salesman 
problem is an example of an integer program for which 
very large instances have recently been solved. The min- 
imum charge formulation of the centrosymmetric phase 
problem most closely resembles a particular class of MIP 
known as the multidimensional 0-1 knapsack problem, a 
A/'P-hard combinatorial optimization problem |l2] , pl| ]. 

Even for a modest-sized problem like i(AlPdMn), the 
search tree for 500 measured reflections (unknown signs) 
has 2 500 leaves, so search by enumeration is impossible. 



Several algorithms have been developed in an attempt 
to solve problems of the kind specified by Eqns. (||)- 
(H). One such algorithm, known as "branch and bound" 
(BnB), has been used to solve some of the largest trav- 
eling salesman problems to date |^|. The basic idea of 
BnB is to use information about bounds on the objective 
function to reduce the search space. 




FIG. 1. Reconstructed electron density in the 5-fold plane 
of the 6D unit cell. Darker regions correspond to higher 
charge density and can be identified with the three types of 
atomic surface proposed by the authors of Ref. 17. 



The feature of integer and mixed integer programming 
problems that makes them so difficult is that the opti- 
mization is over a nonconvex set. If we were to relax the 
integer constraints on the signs and allow s q G [—1, 1], 
we have a related problem called a linear program (LP). 
Linear programs are not A/'T-'-hard so there exist a va- 
riety of polynomial time algorithms to solve them |p3| . 
Though the LP for the phase problem is unphysical, LP 
solutions of integer and mixed integer programs (known 
as "relaxations") serve as upper (lower) bounds on the 
unrelaxed integer programs whose goal is maximization 
(minimization). We may use this information to reduce 
the search space in the following way. Suppose we are at 
depth d in the tree, meaning we have fixed the first d — 1 
signs at ±1 values. Suppose also we have some informa- 
tion about the upper bound on the charge, whether it be 
via feasible though non-optimal solutions, or through a 
priori knowledge. We then solve two relaxed LP prob- 
lems; one in which Sd — 1 and one for which Sd = — 1, 
the relaxation occurring on all signs with index greater 
than d. If one of the relaxed problems yields a value for 
the charge that is greater than our current best value, we 
know the solution cannot lie in that portion of the tree 
and we may prune that entire branch. 
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FIG. 2. Contour plot of the i(AlPdMn) electron density 
within an icosahedral mirror plane centered on the Mackay 
cluster (a center of icosahedral symmetry). The darker re- 
gions correspond to regions of higher electron density. Pairs 
of atoms at unphysically short separation are manifestations 
of partial site occupation in the average structure. The units 
on both axes are A. 



The BnB search strategy that has proven most effec- 
tive is a variation of the "breadth-first" search, with 
the branching order specified by the Fourier ampli- 
tudes. Thus signs corresponding to large amplitudes are 
branched on first, while the very low amplitude signs re- 
main relaxed throughout the optimization. A basic pa- 
rameter of the search is an imposed upper bound on the 
charge. Beginning at the root, the tree is searched ex- 
haustively by branching on successively smaller ampli- 
tude signs until the LP relaxations exceed the imposed 
upper bound. This limits the depth of search into the 
tree, with deeper searches corresponding to higher im- 
posed bounds. Our experience with i(AlPdMn) repro- 
duces the behavior seen in studies with simulated diffrac- 
tion data, where a particular branch of the search tree 
is explored much deeper than all other branches (given 
our imposed bound on the charge). Were one to increase 
the imposed bound, eventually this distinguished branch 
would be the first to reach the top of the tree and yield 
the solution of the optimization problem. This is hardly 
necessary, however, since the important characteristics 
of the solution (large amplitude signs) have already been 
determined when the distinguished branch first appears. 

The efficiency of the algorithm is improved by reduc- 
ing the number of constraints in Eqn. (0). Since each 
inequality is associated with a point of the grid Q, we 
wish to choose the smallest Q consistent with the Fourier 
sampling theorem which applies to our particular space 
group. If G is the corresponding point group, Eqn. (|^) 
generalizes to 
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where $ is the phase function of the space group |25| and 
0(G q ) is the order of the subgroup of G that fixes any 
particular element of the orbit, q 9 . The columns of F rq 
thus correspond to symmetry orbits of reflections. Wc 
make F square and invertible by choosing the number of 
grid points to equal the number of reflection orbits in the 
data set (plus one additional point for q = 0). The lo- 
cations of the grid points within the fundamental region 
of the space group are determined by making F well- 
conditioned. We use a simple procedure which starts 
with random locations and performs conjugate gradient 
displacements that maximize logdetF*F p4j . 

We now present our results for the structure of 
i(AlPdMn), with space group F53^. Due to space con- 
straints we give a detailed report elsewhere [^(|. The 
minimum charge phase determination was performed 
with 499 of the 503 symmetry incquivalent reflections 
in the data set of Boudard et al. [[l7| . The deepest BnB 
search went to depth 84 and required 2730 LP evalua- 
tions. A distinguished branch already appeared in shal- 
low searches with depth as small as 10. As a 6D structure 
i(AlPdMn) is relatively simple and the phase determina- 
tion effort is not so much that of finding the solution but 
more the process of refining (extending) the solution to 
the weaker reflections. 

To minimize on interpretation we only give the re- 
constructed electron density, that is, we do not attempt 
to identify atomic surfaces with idealized polyhedra and 
specific chemical occupation domains. Also, we note that 
since our electron density is derived from Bragg reflection 
data, it represents the average structure and may there- 
fore possess unphysically short interatomic distances, or 
other manifestations of disorder. Strong indications of 
partial/mixed occupational disorder were already found 
prior to phase determination in the process of normaliz- 
ing the experimental data [^7f . 
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The main features in the electron density are seen in 
a "5-fold plane" , Figure |l| — a two dimensional section 
that is commensurate with the 6D lattice and invariant 
with respect to a 5-fold rotation. The compact concen- 
trations of charge in this plane are consistent with the 
three atomic surface model derived by Boudard et al. |i~7| ] 
using Patterson analysis. There is also some evidence in 
our reconstructed "surfaces" of features that go beyond 
the naive tiling models of quasicrystalline structure. For 
example, deviations from a straight "cigar" shape indi- 
cate context dependent shifts from ideal tile vertices. 

By far the most conspicuous feature of the recon- 
structed atomic surfaces is their rounded character. We 
believe this correctly reflects the strong damping of inten- 
sities with the perpendicular component of the 6D wave 
vector and is not an artifact of Fourier truncation. A 
definitive test will require a larger x-ray data set. The 
rounding of the atomic surfaces, particularly in the three 
dimensions perpendicular to physical space, is expected 
in the random tiling scenario |pgj| , where the average 
structure is effectively a superposition of many tile rear- 
rangements. Boundaries of the averaged atomic surfaces 
(in the perpendicular dimensions) would not be sharp 
leading to numerous partially occupied sites and appar- 
ent unphysically short interatomic separations. Our re- 
constructed electron density in a mirror plane through 
the origin in physical space, Figs. || and 0, illustrates 
these effects. Well defined atomic positions occur in the 
second shell and correspond to the Mackay cluster. The 
sites of the first shell, at vertices of a small dodecahedron, 
are only partially occupied to avoid short interatomic 
separations. 

We have described a new algorithm for direct method 
phase determination based on optimization ideas from 
operations research. Our application of the method to 
i(AlPdMn) shows the importance of statistical disorder 
in the structure. We are currently studying the complex- 
ity of phase retrieval in general, including the effects of 
diffraction data resolution. In addition, we are applying 
the method to other instances of phase retrieval, includ- 
ing two-dimensional images and structure determination 
of small biological molecules. 
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